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ABSTRACT 

We study the effects of the initial conditions of turbulent molecular clouds on the ionization structure 
in newly formed H n regions, using three-dimensional, photon-conserving radiative transfer in a pre- 
computed density field from three-dimensional compressible turbulence. Our results show that the 
initial density structure of the gas cloud can play an important role in the resulting structure of the 
H II region. The propagation of the ionization fronts, the shape of the resulting H II region, and the 
total mass ionized depend on the properties of the turbulent density field. Cuts through the ionized 
regions generally show "butterfly" shapes rather than spherical ones, while emission measure maps 
are more spherical if the turbulence is driven on scales small compared to the size of the H II region. 
The ionization structure can be described by an effective clumping factor £ = (n) ■ (n 2 ) / (n) , where n 
is number density of the gas. The larger the value of f, the less mass is ionized, and the more irregular 
the H II region shapes. Because we do not follow dynamics, our results apply only to the early stage 
of ionization when the speed of the ionization fronts remains much larger than the sound speed of the 
ionized gas, or Alfven speed in magnetized clouds if it is larger, so that the dynamical effects can be 
negligible. 

Subject headings: H II regions: ionization structure - H n regions: morphology - radiative transfer - 
ISM: clouds - ISM: hydrodynamics - magneto-hydrodynamics (MHD) - turbulence 



1. INTRODUCTION 

H II regions are photoionized regions surrounding 
young OB stars. The ionization structures and physical 
properties of these regions are important to understand- 
ing of the formation of massive stars and their feedback 
to the environment. Over the past several decades, H II 
regions have been extensively observed in various size and 
shapes (e.g. see reviews bvlWo od fc C hurchwelllll989t 
iCarav fc Lizanolll999t iChurchwellll 19991 120021) "Tjenend- 
ing on their sizes, they are usually classified as ultracom- 
pact (linear size below 0.1 pc), compact (0.1-1 pc), or 
extended (several to tens of parsecs). Ultracompact and 
compact H II regions are ionized by massive stars still cm- 
bedded in their natal molecular clouds, while extended 
H II regions are thought to be in their mature states, pow- 
ered by a combination of stellar radiation, stellar wind 
and supernova explosions in the associated OB star clus- 
ters (Yang et al. 1996). Despite the diversity of sizes, 
H II regions generally display common features such as 
inhomogeneities and irregular shapes. One interesting 
and still open question on the origin of these complex 
structures is: are they from the initial conditions of the 
gas before ionization, or formed by dynamical processes 
during H II region evolution? 

Early analytical and numerical work focused on the 
formation of H II regions and the propagation of ion- 
ization fronts (I- fronts), by modeling the ionization of 
a massive star suddenly turned on in a uniform gas 
(lStromgrenlll939t lKahrill954l lAxford l l 19611 iGoldsworthvl 
II 9611 iTenorio-Taglelll 9761 lElmegreenlri97fiTl . It was soon 
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realized, however, that the interstellar medium is far 
from homogeneous. Simple condensations or smooth 
variations of density were then included in some model 
calcu l ations llFloweill969l:lM7rslJl97fllKirkDatrickll97ol 
[19721 iPeauignot. Stasinska fc Aldrovandil 119781 iKopperl 
1979Hlckdll979llTenorio-Tagle fc Bediidll98U see rYorkel 
1986 for a review). However, in the early studies, 
radiative transfer effects were either ignored in cal- 
culations or limited to one or two-dimensio nal treat- 
ment s, and suffered from low resol u tions ijg andford 
Il97l iKlein. Stein fc Kalkofenl Il973 iKoppenl lL97St 
lTenorio-Taglelll979t iSandford. Whitaker fc Kleinlll982j) , 
which limited their applications to observations. 

Recently, several authors have reported two- 
dimensional simulations of the dynamical evolution 
of H II r egions and their applic at ions to observa- 
tions. iCarcia-Segura fc Franco! l|1996l) presented 
two-dimensional gasdynamical simulations of the evolu- 
tion of H II regions in constant and power-law density 
profiles, and argued that a thin-shell instability of the 
ionization fronts could produce the irreg ular and bright 
edges of the H II region. iWilliama 1)1999!) also performed 
two-dimensional gasdynamical simulations and showed 
that shadowing instab ility could also lead to fo r mation 
of dense clumps. iFrever. Hensler. fc Yorkei ( 2003) 
presented two-dimensional radiation gasdynamical sim- 
ulations of the interaction between an isolated massive 
star and its homogeneous and quiescent ambient ISM, 
and showed that the redistribution of mass by the action 
of the stellar wind shell could form the "finger-like" 
shapes. 

tabel@aHe»pffi£du since star forming regions display 
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a wide range of kinematic properties and den- 
sity distributions, the gas may well be turbulent, 
clump y, filamenta r y, and coupled with the magnetic 
field iElmegreenl Il993t iDvson. William s fc Redman! 
T993 IRedrnan et a.ll fe* IMac Low! Il999t 
Crutcherl 119991: iWilliams. Dvson & Hartauistl 
2000; iKlessen. Heitsch. fc Mac Low! 120001 

Ostriker. Stone, fc GammieT " 120011 see 
Mac Low fe Klessenl 120041 for a recent review). These 
uniform models may not capture the realistic density 
field of the ambient medium. The initial inhomogencitics 
and turbulent structures of the gas clouds may affect 
the propagation of the ionization fronts and contribute 
to the shaping of the resulting H 11 regions. 

To test this, a reliable, three-dimensional treatment of 
radiative processes is essential, and high resolution simu- 
lations of radiative transfer fully coupled with gasdynam- 
ics in turbulent medium are highly desirable. However, 
incorporating radiative transfer in gasdynamics simula- 
tion presents a serious numerical challenge. The difficul- 
ties are due to the non-locality of the radiation physics, 
which hampers efficient parallclization of the simula- 
tions and severely limits the size of problem that can 
be tackled. So far most calculations have not been 
able to include transfer of ionizing radiation into three- 
dimensional gasdynamics simulations, so, as a step for- 
ward to investigate the ionization of turbulent clouds, 
we apply radiative transfer to a static, pre-computed, 
three-dimensional, turbulent density field. It is a sim- 
plified model, but so long as the speed of the I-front 
is much larger than the sound speed of the ionized 
gas, t he expansion of the gas can be neglected l)Spitzerl 
1978), and the density field can be treated as static. 
In cases where a magnetic I-fron t is considered (e.g. 
Willi ams. Dvson fc Hartauistll2000]) , the Alfven speed in 
the mostly neutral gas ahead of the I-front should re- 
place the ionized sound speed in the above criterion if 
it is larger. However, usually the ionized sound speed 
dominates. 

As we will show, this simple approach produces a vari- 
ety of interesting morphological structures that may pro- 
vide us some understanding of the physical processes in 
H 11 regions even before more sophisticated calculations 
can be performed. However, our results only apply to the 
initial photoionization of static turbulence before gas- 
dynamics becomes important, so future comprehensive 
simulations with a full treatment that couples radiative 
transfer and gasdynamics will be necessary to confirm 
and extend them. 

We present high resolution computations of the initial 
ionization of compressible turbulence by a point ionizing 
source, such as a young star. We apply INONE, a t hree- 
dimen sional radiative transfer code developed by lAbell 
(2000), to a compressible turbulent medium from dif- 
ferent hydrodynamical (HD) and mag netohydrodynam- 
ical (MH D) models usi ng ZEUS-3D llStone fc Normanl 
1992© , as described in IMac Low! l| 19991) . INONE inte- 
grates the jump condition of I-fronts, and uses a photon- 
conserving method that is independent of numerical res- 
olution and ensures correct propagatio n speeds for the 
I-fronts (jAbel. Norman fc Madaulll999l) . In 32 we briefly 
describe the methods of computation, and the codes and 
models used in this work; we present the results of propa- 
gation of I-fronts, emission measure maps, and ionization 
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structures in fJ3 and in 0we give an analytical discus- 
sion of the key parameters that determine the ionization 
structure, as well as the validity of the assumptions. 

2. COMPUTATIONS 
2.1. Compressible Turbulent Molecular Clouds 

Star forming molecular clouds have been observed to 
be highly clumpy and have broad emission line w idths, 
which suggests supersonic motions in the clouds (Blitz 
1991. Th ese supersonic motions seem not to be ordered. 
Crutcherl l)1999fl presented Zeeman observations of mag- 
netic fields in molecular clouds that showed that veloc- 
ities in the observed clouds are typically Alfvenic. So 
molecular clouds are generally self-gravitating, clumpy, 
magnetized, turbulent, compressible fluids. Although 
numerical simulations of transient, compressible turbu- 
lent molecular clouds have been carried out for several 
years, we are far from drawing a comprehensive picture 
o f the comp l icated structure of the clouds. 

IMac Low! l)1999|) conducted direct numerical compu- 
tations of uniform, randomly driven turbulence with 
128 3 resolution using the ZEUS-3D MHD code , which 
successfully simulated the density distribution of the 
turbulent medium. ZEUS- 3D is a well-tested, Eule- 
rian, finite-difference c ode l|Stone fc Normanl ^^HEi. 
Clark efc Normanl Il994'). It uses second-order Ivan Leerl 
l)1977[) advection and resolves shocks using von Neumann 
artificial viscosity. It als o includes magnetic fi el ds in 
the M HD approximation l|Hawlev fc Stondll99fllK iFallel 
(2002) pointed out some problems with ZEUS: rarefac- 
tion waves often break up into a series of jumps, and 
adiabatic MHD shocks sometimes show errors. However, 
as we are using isothermal MHD, and are primarily con- 
cerned with the density field, which is barely affected by 
rarefaction shocks, our models remain valid. 

For our computations, we assemble 384 3 de nsity fields 
by re peating the periodic 128 3 models of IMac Low! 
(1999). Table 1 describes the models we use. The model 
names begin with either "H" for hydrodynamic or "M" 
for MHD, then have a letter from "A" to "E" specify- 
ing the level of energy input E{ n , then a number giving 
the dimcnsionless wave-number k chosen for driving, and 
then, for the MHD models, another number indicating 
the initial field strength specified by the ratio of Alfven 
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speed to the sound speed. In the last column, C is the 
gas clumping factor defined as C = (n 2 ) / (n) 2 . Note that 
the values of C listed here are averaged over the whole 
simulated volume. 

2.2. Ionization Front Tracking 

In general, a study of the evolution of ionization zones 
around a ionizing source re quires a full so lution to the 
radiative transfer equation (Kirchhofi 1860): 



1% 

c dt 



(1) 



where n is a unit vector along the direction of the radi- 
ation, I„ is the monochromatic specific intensity of the 
radiation field, and r\ v and Xv are emission and absorp- 
tion coefficients, respectively. 

However, a direct solution of equation Jy| is usu- 
ally impractical becaus e of i ts high dimensionality. 
lAbel. Norman fc Madaul l|1999(l made the calculations 
feasible by developing a ray tracing algorithm for radial 
radiative transfer around point sources, reducing the di- 
mensionality of the transfer equation to a level where 
ionization can be computed on a Cartesian grid. This 
algorithm conserves energy explicitly and thus gives the 
right speed of I- fronts, but ignores the diffuse field pro- 
ducedby scattered radiation. 

Kahn (1954) defined the nomenclature of R ("rar- 
efied" ) and D ( "dense" ) fronts according to the speed of 
the I-fronts with respect to the sound speed of the ion- 
ized gas: R-type fronts move supersonically while D-type 
fronts move subsonically. If a source of ionizing radiation 
suddenly switches on, the I-front is initially weak R-type. 
Once its velocity drops to about twice the s ound speed in 
the io nized gas, the I-front becomes D-type ijDvson et al.l 
2002). The propagation of an R-type I-front in a static 



ctB47rr 2 n p (r)n (r) dr . (2) 



medium is given as: 
^Rfn B (R I )^-=F- 



where R\ is the radius of the I-front, F is the ionizing 
photon flux, 7Jh, Hp and n c are the number density of the 
neutrals, ions, and electrons, respectively. We assume 
a constant case B recombination rate, ob = 2 x 10~ 13 
cm 3 s _1 for T ~ 10 4 K. Integrating equation @ along 
the rays with the ray-tracing technique mentioned above 
gives the time at which the ioni zation front arrives at a 
given cell, as is done in INONE (Abel 2000). Retrieving 
the arrival time in a 3D array allows one to investigate 
the time dependent morphology of the I-fronts. 

2.3. Scaling 

In this paper, we use INONE to compute photoioniza- 
tion of different turbulence models. For simplicity, we 
consider hydrogen gas only. We expect the inclusion of 
helium to have only a minor effect on the results pre- 
sented here. The simulations are scale free and depend 
only on the ionizing flux F and density n of the medium, 
so one can derive the scaling of physical parameters such 
as ionized mass M, ionized volume V, and ionization 
time t in terms of F and n. As the analytical deriva- 
tion in i)4. II demonstrates, the ionized mass at late times 
scales as M oc F/n; the ionized volume at late times is 
close to the volume of the Stromgren sphere, and scales as 



V oc F/n 2 ; and the characteristic timescale is the recom- 
bination time, which scales as t oc 1/n. As an example, 
if we cho ose Fn = 1.2 x 10 49 s" 1 , which is typical for an 
06 star ()Panagial lT973). and a cloud with a size of 0.5 



pc and average density of no = 5 x 10 3 cm" 

r.3 



3 , we then 

have M = 10M Q , V = 0.08 pc 3 , and t Q = 10 11 s (these 
will be given in $3] as the results of the simulations). For 
given F and n, one can rescale our results as follows: 
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In the simulations, we use the pre-computed gas den- 
sity distribution with ZEUS-3D for the HD and MHD 
models listed in Table 1 as input into INONE. We will 
discuss the validity of this treatment later in A4.2I Three 
simulations are conducted for each model, by varying the 
position of the star to the maximum density, minimum 
density, or center of the density field (random density), so 
for our example a total of 30 simulations were performed 
on SGI Origin2000 computers. Each computation took 
several days on a single processor. 

3. RESULTS 
3.1. Shapes of the Ionized Regions 

Figures ^ and [21 show the propagation of the I-fronts 
into the turbulent medium surrounding a point source. 
The background image is the density distribution of the 
gas simulated by HD and MHD models, respectively. 
The contours are the arrival times of the I-fronts. In Fig- 
ure ^ the density fields are taken from HD models with 
different energy input and driving wave number, while 
in Figure [3 the density fields are from MHD models 
with different driving wave number and magnetic field 
strength, as listed in Table 1. In both figures, the left 
panels are for the case that the point source (star) is 
placed at the maximum density, while in the right pan- 
els the star is placed at the minimum density. In each 
individual panel, the contours give the position of the I- 
fronts in our example scaling from 0.1 to 100 years, with 
the interval increasing evenly by a factor of 10. The mod- 
els chosen in each figure have different average clumping 
factors (see Table 1). 

At early times, the recombination term in equation J5J 
can be neglected, and the I-front velocity is determined 
by the local density structure of the gas: the front travels 
faster into voids and slower into denser filaments. The 
shape of the I-front depends on the size and contrast of 
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Fig. 1. — Propagation of an R-type I-front in a 384 3 density field 
simulated with different HD models. The ionizing source is located 
at maximum {left), and minimum density {right), respectively. The 
contours give the position of the I-front from 0.1 yr to 100 yr, with 
the interval increasing evenly by a factor of 10. The size of box 
is 0.5 pc, (no) = 5 X 10 3 cm -3 . Note the image and contour are 
slices through middle point of the cube. 




Fig. 2. — Propagation of the I-front in a 384 3 density field simu- 
lated with different MHD models. The ionizing source is located at 
maximum {left), and minimum density {right), respectively. The 
contours give the position of the I-front from 0.1 to 100 years, with 
the interval increasing evenly by a factor of 10. The size of box 
is 0.5 pc, (no) = 5 X 10 3 cm -3 . Note the image and contour are 
slices through middle point of the cube. 



the voids and clumps. Longer turbulent driving wave- 
lengths (smaller driving wave numbers) create structure 
in the density field on larger scales and produce more 
asymmetric H n regions with characteristic 'butterfly' 
shapes. At later times, recombination becomes signif- 
icant, increasing the dependence on local density and 
making the resulting asymmetry more pronounced. 

The initial I-front velocity is sensitive to the placement 
of the ionizing point source, as can be seen from a com- 
parison of the innermost contours in the left-hand and 
right-hand columns of Figures and [3 However, this ef- 
fect becomes less important at later times, and the final 
H II regions of the same model have very similar sizes, 
regardless the position of the ionizing source, as can be 
seen in Figures [BJ 

In Figures [2 and □ the images are sorted in order of 
increasing clumping factor C; e.g. in FigurenC(HA8) = 
1.5, C(HC2) = 5.88, C(HE2) = 7.24, and similar trend 
in Figure [5] Comparison of the images clearly shows 
a strong link between the size of the clumping factor 
and the morphology of the H n region - the larger the 
clumping factor, the more asymmetric the H n region 
becomes. 

3.2. Emission Measure 

For comparison with observations, we can map the 
emission measure in our models, 

EM = [ n 2 e dl (6) 



where I is the size of the region along the line-of-sight. 
Typically more compact H n regions have higher emission 
measures. For example, an extended H n region with size 
of 10-100 pc and density of 10 cm~ 3 has EM in the range 
of 10 3-4 pc cm -6 , while an ultracompact H n region with 
size < 0.1 pc and density > 10 4 cm -3 is brighter, with 
EM > 10 7 pc cm- 6 . 

We calculate emission measure maps from the simu- 
lations by integrating density n 2 H of the ionized regions 
through the path length in the third direction. Figures 
|2| and 0] show emission measure maps in the x-y plane of 
HD and MHD models, respectively. The arrangement of 
the models in both figures are the same as in Figures ^ 
and We can see that the range of the EM is around 
10 6-8 pc cm -6 for our adopted parameters, which agrees, 
as expected, with the observed range of compact H II re- 
gions. 

The shape and structure of the EM map depends on 
the clumping factor: the smaller the clumping factor is, 
the more spherical the map appears. Filamentary struc- 
tures and contrast within the H n region also increase 
with C. There are some differences between the HD and 
MHD models with C > 1: HD models tend to have more 
irregular shapes, and bigger voids and clumps. Larger 
driving wavelengths of the turbulence also produce larger 
voids and clumps. These results generally suggest that 
the turbulent structures observed in H II regions may be 
driven at large scale, and that the initial gas clouds may 
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Fig. 3. — Emission measure of the ionized regions of HD models 
at t ~ 100 yr, assuming the scales used in the text as an example. 
The ionizing source is located at maximum {left), and minimum 
density {right), respectively. 




Fig. 4. — Emission measure of the ionized regions of MHD models 
at t ~ 100 yr assuming the scales used in the text as an example. 
The ionizing source is located at maximum {left), and minimum 
density {right), respectively. 



have a large clumping factor. However, the influence of 
dynamics will have to be computed to put any such con- 
clusion on a firm basis. 

3.3. Ionization Structure 

Figures and show the time-dependent ionized mass 
and ionized volume, respectively. Note that the units, 
Mo, Voi and to are the same as in equations ©, JU, and 
Different panels have the ionizing star at different 
positions, while within each panel, different curves indi- 
cate different models. The ionized volume is calculated 
by integrating the ionized cells at time t, and the ion- 
ized mass is the integral of the mass encompassed in the 
ionized volume. The curves become saturated around 
t ~ fOO years. It appears that the mass of the ionized 
gas in our example is clos e to that of a compact H II 
region l|Franco et al.1 1200(1) . For the same input energy 
strength and driving wave number, no big difference is 
seen between the MHD and the HD models. However, 
all the curves follow a sorted order of the clumping fac- 
tor C, with HA8 and MA8I having the smallest C, while 
HE2 has the largest C. 

3.4. Compared With Two-Phase Clumpy Model 

Here we compare our results with an analytic model of 
a two-phase clumpy medium. Consider a medium with 
optically thin clumps of number density n c \ in a uniform 
de nsity field of n umber density rib g , which was studied 
bv lKoppenl l )I979j) . If we take a filling factor e for clumps, 
and the density contrast £ = ny, g /n c i, then the average 



density is: 

n = en c \ + (I - e)n hg = n c \ [e + (1 — e)£] , (7) 
The ionization structure of this medium is: 



1 dt 



Ri 



4irr oiBn ,edr 



Ri 



47rr 2 a B nL(l - e)dr , (8) 



where as is the recombination-rate coefficient as in equa- 
tion J2J). Define a two-phase clumping factor C2, 



C 2 = 



e£ + (l-6)£ 2 



(9) 



and define dN = AnR^ndRi, where N is the total number 
of ions. Equation JBJ can then be written as 



dN_ 
~dt 



F - asC^nbgiV , 



(10) 



from which we get an analytical solution of the number 
of the ions N, and the ionized mass M^- 



M 2 = m p -N 



F 



-(1-e 



)■ (11) 



OB^nbg 

Figure [7] shows the ionized mass from the two-phase 
model, and a comparison between the analytical results 
and the numerical ones. The clumping factors C2 were 
chosen to be close to the values measured in the tur- 
bulence simulations, as listed in table 1. From top to 
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Fig. 5. — Ionized mass over time for different models, with ion- 
izing star at (top:) maximum density, (middle:) a random posi- 
tion, and (bottom:) minimum density. Mq and to are denned in 
H2.3I The vertical d ashe d line indicates the minimum time when 
vn = 2c s (see Figure [TuV 
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Fig. 6. — Ionized volume over time for different models, with 
ionizing star at (top:) maximum density, (middle:) a random po- 
sition, and (bottom:) minimum density. Vq and to are defined in 
The vertical d ashe d line indicates the minimum time when 
vn = 2c s (see Figure [TuV 



bottom in the plot, we have 1.0 (uniform gas), 1.5, 3.0, 
4.0, 5.0, 6.0 and 7.3 respectively. The background den- 
sity ribg is taken as the mean density n = 5x 10 3 cm~ 3 
of the scaled simulations. Though similar in shape, the 
results of the two-phase model are distinct from those of 
the turbulent models; the ionized mass for a given Gi is 
less than in the corresponding turbulent model, and the 
difference increases as the medium becomes more clumpy. 

4. DISCUSSION AND SUMMARY 

We have presented simulations of different models, and 
found that a two-phase clumpy model can not fully de- 
scribe the ionization structure of the turbulence mod- 
els. However, some questions still remain: what deter- 
mines the ionization structure? When are our assump- 
tions valid? 





C1.0 0.2 Z:i D.6 0.E 

Fig. 7. — Left: ionized mass of two-phase clumpy model over 
time for different two-phase clumping factors C2 taken after the 
values of the turbulence models in Table 1. Right: comparison of 
numerical ionized mass M n (the ionizing star is located at random 
position) to analytical ionized mass from the two- phase models M2- 
Note that the dotted line in this plot is not a fitting curve, but just 
a diagonal line to guide the eyes. 



4.1. What Determines the Ionization Structure? 

The results shown above suggest that the ionization 
structure is sensitive to the input energy strength and 
driving wave number of the turbulence. But qualita- 
tively, what is the key factor that determines the ioniza- 
tion structure? 

Assume the density field is constant, and the gas is 
fully ionized within the ionization front, then both the 
ion density n p (r) and the electron density n c (r) at radius 
r are equal to the original hydrogen density n(r), n p (r) = 



n c (r) — n(r). Let N(t) be the total number of the ions 
at time t, N(t) = J Q 1 n(r) dV, then equation can be 



rewritten as: 

dNjt) 
dt 



= F-a B I n(r) 2 dV 




F - a B ■ (n) 



(n) 2 



N(t) 



(12) 



We can now see that the ionized mass depends ex- 
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Fig. 8. — Effective clumping factor f versus time as an indicator 
of ionization structure. It varies with different models, and different 
position of the ionizing star at (top:) maximum density, (middle:) 
a random position, and (bottom:) minimum density. 
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Fig. 9. — Comparison of numerical ionized mass M n to analyti- 
cal ionized mass M a for different models, and different position of 
the ionizing star at (top:) maximum density, (middle:) a random 
position, and (bottom:) minimum density. Note that the dotted 
lines in the plots are not fitting curves, but just diagonal lines to 
guide the eyes. 



plicitly on the average gas density (n), and the local 
clumping factor (C = {ri 2 H ) / (tih) 2 ) within the I-front. 
Different driving strength and wave number of the tur- 
bulence yield different size and contrast of clumps and 
voids, which is represented by C and (n), thus yield- 
ing different ionization structures. Since R\ changes 
with time, we further define an effective clumping fac- 
tor ((t) — (n) ■ (n 2 )/(n) 2 , so that equation (jT2|l can be 
rewritten as: 

*?Q=F-a B .{(t).N(t). (13) 

From this we can see that C describes the density field, 
and determines the ionized mass. Figure [5] shows £ as 
a function of time for our models. Comparing Figure |H1 
with Figure \5\ we find that the order of the models in 
Figure |S] is exactly the inverse order in Figure |SJ The 
larger the C, the less mass is ionized, as predicted from 
equation 1131 

The dependence of the effective clumping factor 
on the density field is nontrivial because the averaging 
volume is itself a time-dependent function of the density 
field. However, £ tends to saturate once the radius of the 
I-front is larger than the turbulence driving scale. This 
also corresponds to the saturation of the photoioniza- 
tion. The total ionized mass can then be approximated 
well with equation^! Note also that in these figures, dif- 
ferent positions of the ionizing star does not make much 



difference, because £ is a measure of the density field, 
and the ionized mass depends only on £ and the ioniz- 
ing flux, no matter where the source is. In turbulence 
simulations, resolution might affect the value of £, but it 
would not affect our result as a function of £. 

Figure [5] shows a comparison of numerical and ana- 
lytical values of ionized mass for the different models. 
The small variance between numerical simulations and 
analytical calculations suggests that the calculations are 
self-consistent, and that the numerical resolution of the 
calculation is good enough to resolve the structure. Since 
in reality the density fields are more complicated than 
the simplified model used in the analytical calculations, 
there is a systematic difference between the numerical 
results and the analytical ones: in all these models, the 
smaller the clumping factor, the smaller the difference. 

4.2. When Are the Assumptions Valid? 

In the above simulations, the density field of the gas 
is simulated beforehand and then input to INONE, so 
gasdynamics is not computed simultaneously during the 
course of radiative transfer. However, when the speed 
of the I-fronts slows to roughly twice the sound speed 
of the ionized gas, the expansion of the gas will domi- 
nate the evolu tion of the reg ion, and gasdynamics must 
be considered (Spitzer 1978). Figure ITUI shows the time 
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Fig. 10. — Time evolution of vn/c s , ratio of speed of I-fronts Vjf 
to sound speed c s , for different models. When ti;f decreases to twice 
of the sound speed, expansion of the ionized gas is significant. The 
vertical dashed line indicates the minimum time at when V\{ = 2c s ■ 



evolution of the ratio of speed of I-fronts Vn to sound 
speed c s , Va/c s , for different models. The sound speed 
is c s = 9.1 km s" 1 for T = 10 4 K. The velocity of the 
I-fronts, Vif = dRi/dt, is hard to derive because it is 
three dimensional, and even at the same photon arrival 
time, it may vary dramatically in different cells if the 
local clumping factors are different, so what is shown 
here is an average value at Ri, the radius of the I-front 
at time t. However, one should keep in mind that in 
some regions vit may be much smaller than the aver- 
age value. For example, an obliquely propagating ion- 
ization front moves slow ly and transitio ns from R-type 
to D-type more quickly jWilliamj ll999l). Fully dynam- 
ical simulations are necessary to assess the impact of 
these effects. From Figure El we can see that the speed 
of I-fronts decreases rapidly with time, and at about 
logio(t/to) = —1-4 (which corresponds to t ~ 100 years 
in our example), it drops to twice of the sound speed, 
so after this time the simulations are likely unreliable. 
From Figure |3] and we can see that this is the time 
when the curves of ionized mass and ionized volume be- 
come saturated, so the details of the ionization structure 
d erived up to that point are like l y relia ble . 

Williams . Dvson fc HartauistJ l)2000|) presented one- 
dimensional simulations on jump conditions of magnetic 
I-fronts and showed that fast R-type I-fronts could de- 
crease quickly to slow D-type I-fronts in oblique magnetic 
fields. We note, however, that as long as the propaga- 



tion speed of the I-fronts is much larger than the Alfven 
speed of the medium, this effect remains unimportant. 
In molecular clouds, the typical Alfven speed is about 
1-5 km s -1 , which is no larger than the sound speed in 
ionized gas. Table 1 shows that the ratio of Alfven speed 
to the sound speed in the MHD models is appropriate 
for cold, molecular clouds. iCrutcheil l)1999fl showed an 
empirical power-law relationship between magnetic field 
strength and density above 10 3-4 cm -3 . According to 
this scaling law, our density field ((n) = 5000 cm~ 3 ) 
might have B ~ 50 /iG, which gives an Alfven velocity 
va = (B 2 /(47Tyo)) 1 / 2 ~ 5 km s _1 . The ionized sound 
speed c s = 9.1 km s _1 , larger than the Alfven speed va, 
which suggests that our critical point logio(t/io) = —1-4 
occurs before the magnetic field becomes important. 

Again, it should be emphasized that our simulations 
are simplified models in which a point ionizing source is 
turned on in a static turbulent density field. We did not 
take into account several effects, such as the mass accre- 
tion by the ionizing star, the diffuse field, and an y time 
evolu tion in the ionizing luminosity of the star l|Yorkel 
1986). These effects should eventually be included in 
more comprehensive future simulations. Our main focus 
here is to investigate how the initial conditions of the 
gas affect the ionization structures of the resulting H n 
region. 

4.3. Summary 

We have presented a study of the effects of turbu- 
lent initial conditions on the ionization structures of 
young H n regions. We performed high resolution, three- 
dimensional radiative transfer computations, as well as 
analytical calculations, of the propagation of I-fronts into 
static density distributions drawn from numerical simu- 
lations of compressible, turbulent molecular clouds. Our 
results show that the initial turbulent structure of the 
gas can play an important role in the early ionization 
structure. The propagation of the I-fronts, the shape 
of the resulting H II regions, and the total mass ionized 
depend on the local strength of the voids and clumps. 
We can characterize this with an effective clumping fac- 
tor £ = (n 2 )/(n). The ionization fronts move quickly 
into voids but slowly into dense filaments. The ionized 
mass changes inversely with £. Cuts through the ionized 
regions generally have "butterfly" shapes, with larger £ 
producing more irregular shapes. Larger scale driving of 
the turbulence produces more filamentary emission mea- 
sure maps. 

We emphasize that our results are based on static den- 
sity fields, so they apply only to the early stages of the 
evolution of an H II region when the speed of the ioniza- 
tion fronts remains much larger than the sound speed 
of the ionized gas, so that the expansion of the gas 
is negligible. (In cases where magnetic field is consid- 
ered, the neutral Alfven speed should be used in addi- 
tion to the ionized sound speed, but it is rarely higher.) 
We regard our results as complementary to t he fin dings 
oflOarcia-Segura fc Francol ( 1996). Williams! (|1999fl . and 
iFrever. Hensler. fc Yorkel (|2003f) . who study the dynam- 
ical evolution of H II regions with two-dimensional simu- 
lations, but only with uniform or smoothly varying initial 
conditions. Future comprehensive simulations with a full 
treatment that couples gas dynamics and radiative trans- 
fer, and with realistic initial conditions, will be necessary 
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to extend these results to later times. 
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